\chapter{\label{chp:intro}Introduction}
\section{\label{sec:small_animal_optical_imaging}Small animal optical imaging}
With the advance of genomic technologies, elaborate animal models with human disease phenotypes enables the study of complex disease patterns and the response of disease to environmental factors and drugs. These advances triggered a major shift of biological research from petri dish and clinical practice to preclinical investigations. The wide variety of endogenous and exogenous contrast mechanisms in small animals allows for visualizing the biochemical, genetic or pharmacological processes \emph{in vivo} and non-invasively, by preclinical imaging modalities. As a central tool in biological research and medical diagnosis, research in preclinical imaging techniques has grown exponentially in the past decade \cite{Lewis2002,Koo2006}.

% small animal modalities
Targeting sensitivity and specificity for retrieving anatomical, functional, and molecular parameters noninvasively in entire animals, a number of clinical imaging techniques have been adapted for animal use, including X-ray tomography, magnetic resonance, nuclear imaging and high frequency ultrasound approaches. Among the numerous imaging techniques, photonic methods provides significant advantages for its high detection sensitivity, low cost, simplicity of implementation, and ability to utilize various contrast mechanisms. Within photonic methods, specific biological processes can be targeted using biomarkers, thanks to the engineering of molecular probes with distinct optical contrast characteristics. These advantages generate great interest in both academia and industry in developing new technologies and applications for optical small animal imaging.

% optical systems..
Latest technologies for small animal optical imaging include the development of targeting agents and the design of imaging systems (as seen in Table \ref{tab:intro-platform}) that allow for observations of protein activity, cellular and subcellular function, gene regulation, physiological responses, and detailed anatomical information. As a noninvasive imaging technique, optical imaging yields longitudinal observations on the same animal as a function of external stimuli, compared to using a large number of animal cohorts generally employed for statistical reliability. Optical imaging can thus serve as a cost efficient tool for diagnostics, monitoring drug delivery and therapeutic interventions in preclinical settings, with clear potential for translation into the clinic \cite{Hielscher2005}.

\begin{table}
\caption{\label{tab:intro-platform}{Summary of the main technical features associated with high-end commercially available whole-body fluorescence imaging systems (CRi=Caliper; VisEN=Perkin Elmer; CW: Continuous Wave; TD: Time Domain).}}
\centering
\def\temptablewidth{\textwidth}
{\rule{\temptablewidth}{1pt}}
\begin{tabular}{cccccc}
 &Company& Domain & Geometry & 3D imaging & Dynamical \\
\hline
KODAK           &Care stream    & CW& Epi-illum.  & ${\times}$ &  ${\surd}$\\
Maestro         & CRi           & CW& Epi-illum.  & ${\times}$ &  ${\surd}$ \\
IVIS Spectrum   &Caliper        & CW& Trans-illum. & Tomography &  ${\times}$ \\
IVIS Kinetic    &Caliper        & CW& Epi-illum.  & ${\times}$ &  ${\surd}$\\
FMT 2500        &VisEn          & CW& Trans-illum. & Tomography &  ${\times}$ \\
iBox            &UV P           & CW& Epi-illum.  & ${\times}$ &  ${\surd}$\\
Pearl Imager    &Li-COR         & CW& Epi-illum.  & ${\times}$ &  ${\surd}$\\
NightOWL II     &Berthold       & CW& Epi-illum.  & ${\times}$ &  ${\surd}$\\
Optix MX3       &ART            & TD& Epi-illum.  & Topography &  ${\times}$
\end{tabular}
{\rule{\temptablewidth}{1pt}}
\end{table}

\section{Imaging optical contrast}
The optical imaging techniques provide fundamentally different information than the other modalities by characterizing the photon-tissue interactions: absorption of light, scattering of light, and emission of fluorescence (Fig. \ref{fig:intro-contrast}). The signal detected by optical imaging techniques is affected by the degree of these interactions, which relates to various state in the physiological status of the imaged sample, and results in the contrasts for optical imaging. In this section, coefficients and characteristics of photon-tissue interactions are introduced and discussed.
\begin{figure}[htbp]
\centering
\includegraphics{img/intro-2}% Here is how to import EPS art
\caption{The photon-tissue interactions in optical imaging: (a) absorption and scattering of light, and (b) fluorescence emission. \label{fig:intro-contrast}}
\end{figure}\\

\noindent {\bf Absorption}

When a photon interacts with a single particle, some of the photon energy may be absorbed and converted into heat. If only the absorption process occurs in a homogeneous media, the transmitted light intensity $I$ after absorption will follow the Beer-lambert's law:
\begin{equation}
I = I_0\exp(-\mu_al),
\end{equation}
where $I_0$ is the incident intensity, $l$ is the optical path length, and $\mu_a$ is known as the absorption coefficient. The absorption of biological tissue is minimal in the near-infrared (NIR) wavelength, where the value of $\mu_a$ is determined mainly by the concentration of oxy-hemoglobin(HbO$_2$), deoxy-hemoglobin(Hb), lipids and water (Fig. \ref{fig:intro-spect}, from \cite{Venugopal2011}). Other substances, such as melanin, as well as several structural tissue components, such as collagen, elastin and lipo-pigments, may have small contribution to it but are generally not considered significant.

\begin{figure}[htbp]
\centering
\includegraphics{img/intro-1}% Here is how to import EPS art
\caption{The absorption spectra of major tissue chromophores. \label{fig:intro-spect}}
\end{figure}

In this spectral window, light can transmit through tissues with thickness of a few centimeters, while for wavelengths below 600nm, high attenuation of light results in a small penetration depth of hundredths of micrometers up to a few millimeters. Hence, at shorter wavelengths, only superficial assessment of tissues is possible, while NIR light allows for imaging preclinical models entirely.

The total absorption coefficient at a specific wavelength $\lambda$ of the tissue can be expressed as a linear combination of the concentration of different chromophores with respect to their extinction coefficients:
\begin{equation}\label{eq:intro-mua}
\mu_a^{\lambda} =
  \varepsilon_{Hb}^{\lambda} [Hb] + \varepsilon_{HbO_2}^{\lambda} [HbO_2] +\varepsilon_{H_2O}^{\lambda} [H_2O] + \varepsilon_{lip}^{\lambda} [lip].
\end{equation}
Notice here the products of the extinction coefficient $\varepsilon$ at wavelength $\lambda$ and the concentration of different chromophores (denoted as $[chromophore]$) are their respective absorption coefficients. This linear relationship can be used to estimate the concentration of chromophores by reconstructing the wavelength-dependent absorption coefficient of the tissue. In this case, multispectral measurements are required and this approach is referred to as the functional optical imaging.\\

\noindent {\bf Scattering}

Photons do not experience pure absorption in most bio-tissues. Scattering of light is another process for photon-tissue interaction that occurs, where a photon changes direction from a straight trajectory upon interacting with particles in the medium through which it passes. This process is considered elastic without energy loss and is the predominant interaction occurring during light propagation in the NIR window. The probability density function of the distance between two consecutive scattering events (step size $l_s$) is given by:
\begin{equation}
\label{eq:prob_scat_length}
P_{l_s}(l_s) = \mu_s\exp(-\mu_s l_s),
\end{equation}
where $\mu_s$ is named the scattering coefficient, which is the inverse of the average step size $l_s$. Like the absorption coefficient, the scattering coefficient is also wavelength dependent and follows the power law \cite{Mourant1997, Torricelli2001}:
\begin{equation}
\mu_s = A\lambda^{-b},
\end{equation}
where $A$ and $b$ represent the scattering amplitude and scattering power, respectively. %need the biological relation to A and b?

When a scattering event occurs, the direction of the scattered photon has to be determined. The Henyey-Greenstein phase function has been commonly used to approximate the probability density function of the direction of photon propagation after an scattering event, by introducing a single parameter scattering anisotropy $g$ \cite{Henyey1941}. $g$ has a value from 0 to 1, where 0 indicates isotropic scattering and 1 indicates forward direction after scattering:
\begin{equation}
\label{eq:henyey}
P_{HG}(g,\theta) = \frac{1-g^2}{4\pi(1+g^2-2g\cos\theta)^{3/2}}.
\end{equation}
$g$ in biological tissue is usually between 0.69 and 0.99 \cite{Cheong1990}. In multiple scattering (diffusion) scenario, it is convenient to consider the total effect of $\mu_s$ and $g$ in the form of the reduced scattering coefficient
\begin{equation}
\label{eq:power_law}
\mu_s' = (1-g)\mu_s.
\end{equation}

\noindent {\bf Fluorescence}

%what is fluorescence
In addition to the absorption and scattering processes, the characteristics of incident photons can be altered by the fluorescence phenomenon as well. When a photon is absorbed by a fluorophore, all the energy possessed by the photon is transferred to the fluorophore in the form of an excited, unstable electronic state of a fluorophore atom. After a certain delay in time, the excited fluorescent molecule will release the absorbed energy and return to the ground state by emitting a photon at a longer wavelength. The difference in wavelength between the excitation and emission spectrum enables visualization of fluorescing objects by filtering out the excitation light. Moreover, the spectra for different fluorophores may vary, thus the light from multiple fluorophores can be distinguished using wavelength filters.

Another way to unmix the signals generated by different fluorophores is to use the lifetime, an inherent parameter of a fluorophore. Lifetime is determined by the time the molecule stays in its excited state before emitting a photon. This delay in time typically follows an exponential rule and reflects the decay of fluorescence:
\begin{equation}
F(t) = F_0\exp(-t/\tau),
\end{equation}
where $t$ is time, $\tau$ is the fluorophore's lifetime, $F_0$ is the initial fluorescence intensity at $t=0$.

%quantify
The fundamental principle for quantifying a fluorochrome is based on the fact that the intensity of the emitted fluorescent light is proportional to the product of the amount of fluorochrome present, its extinction coefficient $\varepsilon$ and its quantum efficiency ($\phi$), when the wavelength and intensity of the illuminating light are constant. The quantum efficiency of a fluorochrome can be highly sensitive to a number of environmental factors, including temperature, ionic strength, pH. etc \cite{Jacques2003, Han2010}. Thus sometimes the effective quantum yield $\eta$, considering both the concentration of the fluorochrome and its quantum efficiency, is employed for quantification. In biological and medical applications, proteins or other cell components can be labeled with extrinsic fluorescent dyes, thus the linear relationship of the effective quantum yield to its emission light can serve as a non-destructive mean of analyzing biological molecules. The functional information it provides by locating and tracking the bio-distribution of the target in tissue is extremely valuable in the drug development process.

\section{Optical tomography}

Optical imaging instrumentation includes a wide variety of technologies, ranging from high-resolution microscopy and optical-coherence tomography for examining tissue at the cellular level, to high-sensitivity noninvasive imaging of deep tissue. Optical methods offer the possibility to image numerous molecular targets with multiple distinct agents similar to immunofluorescence microscopy and fluorescence cytometry. Biophotonics is able to simultaneously distinguish five or more separate imaging probes \emph{in vivo}~\cite{Kobayashi2007, Kosaka2009}. Current \emph{in vivo} fluorescence multiplexing studies are mainly conducted in preclinical settings with Fluorescence Reflectance Imaging (FRI)~\cite{Ntziachristos2005, Sevick-Muraca2008}. However, due to the limited information collected by FRI, the technique is unable to resolve the signal depth and hence quantify it. Moreover, due to the predominance of scattering, planar imaging suffers from low resolution. For these reasons, usage of FRI is limited to superficial observations of subcutaneous flank xenograft models, surgically exposed organs or for intra-operative use, and is not appropriate for studying advanced disease models in internal organs. For such applications, tomography techniques are required.

Optical tomography techniques are based on the characteristics of photons at long wavelength (700-900nm) that are able to probe deeper tissue (a few centimeters). This kind of photons is noninvasive, and particularly effective in small animals where photons can be transmitted throughout the entire body. Similar to x-ray tomography, the optical tomography relies on sequentially exciting multiple points on the tissue boundary and collecting diffused light exiting a few centimeters away from the excitation points at multiple detectors, as illustrated in Fig. \ref{fig:intro-tomo}a. Tomographic techniques utilizing light propagation models thus allow three-dimensional visualization \emph{in vivo} via an algorithm described as follows.
\begin{figure}[htbp]
\centering
\includegraphics{img/intro-3}% Here is how to import EPS art
\caption{(a) Optical tomography source (S) - detector (D) settings for reflectance (top two detectors) and transmission (bottom two detectors) geometries. The dashed arrows show the sources that are sequentially excited at other positions. (b) The weight matrix for a NIR source-detector pair in optical tomography. \label{fig:intro-tomo}}
\end{figure}

\subsection{Forward model}
Like most imaging techniques, the algorithm for tomography can be divided into two stages. First, it requires a forward model that describes the propagation of photons through the medium. If the region of interest is spatially discretized into a set of volumn elements (eg. voxels or tetrahedral elements), this step can be expressed in the form:
\begin{equation}
M(r_s, r_d) = F(r_s, r_d, \alpha(r)),
\end{equation}
where the forward model is denoted as a function $F$ of the source position $r_s$, the detector position $r_d$ and the unknown parameter $\alpha$ of interest, which is a vector variable dependent on the position (voxel) $r$. The predicted measurements is denoted as $M$. For example, in x-ray tomography, the log of measured intensity $I$ divided by the initial intensity $I_0$ is:
\begin{equation}
\ln(I/I_0) = -\sum_{i=0}^N\alpha_jl,
\end{equation}
where N is the number of voxels along the line connecting the source and detector, $\alpha$ is the absorption coefficient at this line as a function of different positions, $l$ is the length of the voxel. The quantity that relates the unknown to the measurement is called weight. In this case the weights along the line are $-l$ and zeros elsewhere, due to the fact that the detected x-rays rarely scatter outside of the line. The weight describes how sensitive a change at this position will result in a difference in the measurement, so the weight map is sometimes referred to as the sensitivity map as well.

The relationship between the absorption coefficient and the log measurement in x-ray is perfectly linear. However, when the forward model is not linear with respect to the unknowns, it is preferable to linearize the forward model to reduce the complexity. Consider a background photon measurement $U(0)$ at a background $\alpha_0$ with a sufficiently small perturbation $\delta \alpha$ so that $\alpha = \alpha_0 + \delta \alpha$, we can obtain the Born expression by a first order Taylor series expansion of $U$ about $U_0$:
\begin{equation}
U = U_0 + \delta U = U_0+\frac{\partial F(\alpha_0)}{\partial \alpha} \delta \alpha.
\end{equation}
The assumed small $\delta \alpha$ is essential in this formula, and this method is known as the perturbation method. Now the changes in $\alpha$ would alter the values of the measurements in a linear way. In this case the first derivative of the forward model with respect to the change in parameter for reconstruction at different position will be the weight to form the sensitivity map.

%The Rytov formalism assumes the heterogeneity perturbs. The total field, ,is therefore given by
%\begin{equation}
%U = \exp(\psi_0+\delta\psi) = U_0\exp(\delta \psi)
%\end{equation}
%The linearization based on small $\delta \psi/\psi$ leads to the Rytov series, on which the first term is given by

Unlike the ``straight line'' propagation for x-ray, the photon propagation in tissue is diffusive due to the predominance of scattering events, resulting in a ``banana shape'' of the sensitivity map (Fig. \ref{fig:intro-tomo}b). To date, analytical models such as the radiative transport equation (RTE) and its approximation, the diffusion equation (DE), have been developed to solve the forward problem. Moreover, numerical solutions using finite difference, finite element methods or Monte Carlo (MC) simulations in heterogeneous media are also employed in calculating the forward problem. These solutions give ways of expressing the weight as functions of the optical properties of the medium and the configuration of the source and detector. The relative merits of these techniques are discussed in the next section.

\subsection{The inverse problem}
The goal of tomographic imaging is to reconstruct the parameters of interest by minimizing the difference between the measured signal and the data predicted by the forward model on the estimate $\alpha$. After the sensitivity function is assigned to each voxel, the unknown parameter can be casted into an inverse problem, that is, the measurements are linear combinations of the sensitivity maps and unknowns:
\begin{equation}\label{eq:linear1}
\left[\begin{array}{c}
 \delta U_1 \\
 \vdots \\
 \delta U_n
 \end{array}\right]
 =
 \left[\begin{array}{ccc}
 \frac{\partial F}{\alpha_{11}} &\dots&\frac{\partial F}{\alpha_{1m}}\\
 \vdots &\ddots&\vdots\\
 \frac{\partial F}{\alpha_{n1}} &\dots&\frac{\partial F}{\alpha_{nm}}\\
 \end{array}\right]
 \times
 \left[\begin{array}{c}
 \delta\alpha(r_1) \\
 \vdots \\
 \delta\alpha(r_m)
 \end{array}\right],
\end{equation}
where $m$ and $n$ enumerate the voxels and the measurements, respectively. The sensitivity matrix in this case, is formed as the Jacobian matrix of the forward model with respect to the unknowns, thus also known as Jacobian in many descriptions.

Now the unknown parameter can be retrieved by inverting the linear system. Popular inversion methods include direct inversion of the matrix and error minimization algorithms. The direct methods may
result in long calculation time for inverting the large matrix (hours). By using the error minimization methods, the unknowns can be solved by fitting the forward model to the measurements very quickly (in minutes or seconds). There are many minimization techniques for solving linear systems, including the algebraic reconstruction techniques (ART) \cite{Intes2004} and the conjugate gradient method \cite{Hielcher1999, Arridge1999}. The objective of these methods are to find the shortest path from the initial guess to the local minimal. In many of the tomographic applications, a least square form of the fit for the predicted detector reading and real measurements is employed as the objective function for minimization.

\subsection{Tomography for different parameters}
The tomography technique has been explored broadly in clinical settings (for breast \cite{Intes2003, leff2008} and brain imaging \cite{Selb2006}) and has been successfully applied to small animal optical imaging to retrieve 3D distribution of various contrasts \cite{Hielscher2005}. Dependent on the difference of the reconstructed quantities, optical tomography methods can be roughly classified as diffuse optical tomography (DOT) and fluorescence molecular tomography (FMT). They can share the same imaging modality and have fundamentally similar principles for reconstruction; however, many differences still exist between the two.

DOT offers unique benefits of direct 3D quantitative recovery of functional tissue properties in thick tissue. Not only can the contrast be brought from the intrinsic optical properties, but also from extrinsic molecules for augmented contrast in absorption and scattering. As the functional parameters are related to the optical properties, this method enables an effective way of estimating the functional status of the bio-tissue in the volume of interest. Moreover, as a non-linear process dominated by scattering, light propagation in tissue is highly dependent on the accurate quantification of the absolute optical properties of the probed specimen, which is essential for rigorous photon propagation modeling in tissue and, hence, important for improving optical reconstruction fidelity \cite{Cheng1999}.

In contrast, FMT aims at obtaining the 3-D fluorescent marker distribution \cite{Mansfield2010}. This makes FMT a suitable method for small animal research where fluorophores can be designed to label the drugs of interest, thus to track their delivery processes. The forward model for fluorescence first includes the propagation of the excitation photons, then the emission and propagation of fluorescent photons from the fluorescing target as a secondary source. When more than one kind of fluorophores are employed in imaging, multiplexing methods should be applied to distinguish the signals from different fluorophores. One method is to collect spectral datasets for efficient unmixing of independent signals~\cite{Mansfield2010}. Such datasets are sequentially acquired, one wavelength at a time, thus lead to relatively long acquisition time, especially for whole-body applications. Alternatively, fluorophore unmixing can be performed based on lifetime contrasts~\cite{Cubeddu2002}. By recording time dependent fluorescence signals, it is possible to distinguish and estimate the fractional contribution of multiple fluorophores with a monochromatic dataset, allowing for fast acquisition protocols. Thus, lifetime based multiplexing is the most promising approach to simultaneously image multiple biomarkers for whole-body applications \emph{in vivo}.

% target: accurate and efficient
\section{\label{} Specific challenges in preclinical tomographic applications}
The implementation of optical imaging techniques in preclinical scenarios faces several unique challenges that need to be addressed, in order to produce robust and quantitative imaging platforms required for widespread acceptance. These challenges include setting up suitable illumination-detection schemes, selecting an accurate forward model, and applying appropriate data in the inverse model.

\subsection{\label{}Forward light propagation models}
%not accurate
To perform quantitative optical tomography, it is essential to employ an accurate mathematical model to describe the propagation of light in biological tissues. However, whole-body small animal imaging is characterized by complex boundaries and complex spatial distribution of small organs that exhibit a wide range of optical properties. In this case, the RTE is not efficient to solve analytically or numerically. The small finite volumes and their associated complex boundary conditions, the broad span of optical properties and potential low-scattering properties of certain organs may limit the validity of DE based models in small animals~\cite{Rasmussen2007}.

%MC accuracy
Alternatively, the Monte Carlo method has been proven to give accurate solusions for the light propagation over a large span of optical properties (or conversely over a large spectral range) and for small interrogated volumes encountered in preclinical scenarios \cite{Wang1995, Ma2007,Rasmussen2007}. A number of studies focusing on Monte Carlo methods for fluorescence signal prediction have demonstrated its accuracy using synthetic and experimental data~\cite{Liebert2008, Swartling2003, Boas2002}. Moreover, when employing time-gated data (will be discussed in Section \ref{sec:intro-domain}), the Monte Carlo approach has major advantage for accurately simulating the transition between minimally scattered photons and diffuse photons (early rising part of the measured temporal point spread functions (TPSFs))~\cite{Mitra1999, Wang1995a, Yaroslavsky1997, Sakami2002, Xu2002, Cai2003, Valim2010}, which DE cannot accurately represent~\cite{Yoo1990, Niedre2010,Yaroslavsky1997}. For example, Valim has shown the time-resolved Monte Carlo simulations agreed better with experimental data than the DE in the comparison of FWHMs for the rising and quasi-CW portion of TPSFs as a response to an impuse source \cite{Valim2010}.

The Monte Carlo method is thus considered the gold standard to model light propagation in either diffusing or non-diffusing media with flexibility to model arbitrary boundaries \cite{Yoo1990}. It has been selected as the most accurate forward model to validate newly developed algorithms (based on the diffusion equation or the radiative transport equation) for fluence/flux prediction, and has been studied for years in spectroscopy and bulk optical property estimation numerically or experimentally \cite{Pan2007,Hayakawa2001, Hayakawa2004, Kienle1996, Palmer2006}. For instance, Hayakawa has experimentally demonstrated that their Monte Carlo based algorithm exhibits the capability of determining $\mu_a$ and $\mu_s'$ for media $0.33\leq(\mu_s'/\mu_a)\leq300$ with errors less than $\sim15\%$. This range of optical properties spanning nearly 3 orders of magnitude distinguishes the Monte Carlo based method from those using radiative transport predictions or the standard diffusion approximation.

% accuracy showcase
A comparison of the Monte Carlo method with the diffusion equation using a publicly available package (Photon Migration Imaging (PMI) toolbox) \cite{PMI} is provided in Fig. \ref{fig:MC_DE} for geometry and optical properties that are mimicking the liver and lung of a small animal. The simulations were performed in a slab geometry to allow for analytical solutions of the diffusion equations. The DE and Monte Carlo forward models in the case of organs such as the liver exhibits a 20\% difference in the light transmitted through the slab, highlighting the inability of DE based algorithms to accurately model light propagation in preclinical scenarios.
\begin{figure}[htbp]
\centering
\includegraphics{img/MC_DE}% Here is how to import EPS art
\caption{Comparison of the light propagation computed by the Monte Carlo method with that by the diffusion equation. The optical properties are at wavelength of 675nm. The size of the inclusion is roughly mimicking the organs in preclinical models. Background optical properties are $\mu_a$ = 0.28cm$^{-1}$, $\mu_s'$ = 14.7cm$^{-1}$, and the optical properties for (a) liver are $\mu_a$ = 1.37cm$^{-1}$, $\mu_s'$ = 6.71cm$^{-1}$, and for lungs are $\mu_a$ = 0.50cm$^{-1}$, $\mu_s'$ = 21.63cm$^{-1}$.}
\label{fig:MC_DE}
\end{figure}

%%drawback;
Nevertheless, Monte-Carlo-based reconstructions are not performed and DE based algorithms are the norm in small animal imaging. Owing to the statistical nature, Monte-Carlo-based techniques suffer from low computational efficiency (hours/days of computation) \cite{Boas2002} because numerous photon trials are required to obtain satisfactory results with adequate signal to noise ratio (SNR). This translates into large computing power thus makes the Monte-Carlo-based methods  computationally expensive. This computational burden, as well as the memory constraints in early computers have hindered the use of Monte Carlo methods as a forward model for tomographic reconstruction, where simulations of numerous source-detector pairs are required and further increase the calculation time in orders of magnitude.

However, there is renewed interest in implementing the Monte Carlo method for tomographic use because of its appealing accuracy in modeling light transport for preclinical studies \cite{Kumar2008a, Raymond2010, Zhang2009b}. For instance, continuous wave weight functions based on Monte Carlo simulations in an adjoint form have been recently applied to tomographic reconstructions in association with unmixing techniques fitting the decaying part of the time-domain data~\cite{Kumar2008a, Raymond2010}.

\subsection{\label{sec:intro-chal-wide}Source and detector configuration}
% point source drawback
Optical tomography performances are critically dependent on the number of spatial measurements collected. To achieve relative high-resolution and hence accurate quantification, dense layout of source detector pairs are required. With the advent of CCD technology, it is straightforward to acquire sufficient detector measurements in parallel. However, current optical tomography systems rely on point excitation sources. To obtain a dense spatial sampling with punctual sources, either dense array of fixed fibers or raster scanning of a laser beam over the surface of the specimen imaged are utilized. The acquisition of measurements for all source positions is then time-multiplexed. For preclinical scenarios, such illumination strategies leads to imaging protocols on the order of hours to image the animal whole-body (neck to tail).

% widefield/ structured
Conversely, FRI, which employs full-field illumination, is a very fast technique but does not collect tomographic data.
However, recent development in micro-mirror technology allows
spatial encoding in widefield illumination to increase the
information content collected. Spatial coding of the widefield
illumination has been successfully employed in fluorescence
microscopy and termed structured illumination. The structured light
methods have been extensively used for profilometric measurements
\cite{Battle1998_PR} and for machine vision \cite{Will1971_AI}.

Recently, such methods have greatly impacted bio-imaging
applications. In microscopy specifically, methods that use
non-homogenous illumination have been shown to improve resolution by
exposing normally inaccessible high-resolution information in the
observed image \cite{Neil1997_OL}. Numerous advances in structured illumination schemes have been achieved, including the introduction of advanced axial and lateral patterns \cite{Fedosseev2005_OLE}, three-dimensional structured illumination patterns and single exposure optical sectioning with color-structured patterns \cite{Krzewina2006_OL}. Transillumination \cite{Pitris2005_OL} and reflectance \cite{Tkaczyk2004_OE} versions of structured illumination have also been introduced. Likewise, the utility of structured illumination for imaging and quantification of optical properties in thick tissues such as biological tissues has been demonstrated \cite{Cuccia2005_OL}.

Modulated Imaging (MI) is based on the same principles as structured illumination microscopy and consists of illuminating the specimen with a widefield, sinusoidal spatially modulated light \cite{Cuccia2009_JBO}. The spatial modulation transfer function of the turbid sample encodes both depth and optical property information, enabling both quantification and sectional imaging of the spatially varying medium's optical properties (Fig. \ref{fig:intro-pattern}, from \cite{Cuccia2009_JBO}). By analyzing the frequency-dependent reflectance, one can quantitatively sample the optical properties of the medium and by varying the spatial frequency of the illumination pattern, one can control the depth sensitivity of detection inside the turbid medium, whereas using planar illumination corresponds to a fixed mean depth of interrogation. The approach has been successfully applied to functional imaging of preclinical models during ischemic stroke \cite{Cuccia2009_JBO}. However, the tomographic capabilities of the technique are dependent on the illumination spatial frequency. The optical sectioning (tomography) is then confined to shallow tissue (a few millimeters) such as thinned area of the rat skull \cite{Abookasis2009_JBO}. To extend the technique to deep seated tissue or organs, true tomographic principles based on inverse formulation should be implemented.
\begin{figure}[htbp]
\centering
\includegraphics{img/intro-4}% Here is how to import EPS art
\caption{Schematic of modulated illumination source (in the $x$ direction only) and the resulting modulated internal fluence rate with the same frequency and phase.\label{fig:intro-pattern}}
\end{figure}

\subsection{\label{sec:intro-domain}Measurement datatype}
The measurement datatype plays an important role in optical tomography. Three optical technology domains exist in the diffuse optical imaging field: continuous wave domain, frequency domain, and time domain, with each generating corresponding type of measurements, as depicted in Fig. \ref{fig:intro-pattern}.\\

\begin{figure}[htbp]
\centering
\includegraphics[width=\textwidth]{img/intro-5}% Here is how to import EPS art
\caption{Illustration of the three optical domains: a) CW, b) FD and c) TD. \label{fig:intro-pattern}}
\end{figure}

\noindent {\bf Continuous wave measurements}

Continuous wave (CW) methods measure the steady state of the light
signals. Instrumentation in CW is the most widespread technique due
to its simplicity to implement and low cost components; however the
CW measurements provide limited information content leading to
non-unique solution in the inverse problem~\cite{Arridge1998}. Thus,
CW techniques are unable to estimates and resolve the tissue
absorption and scattering properties, or fluorescence lifetime \cite{Kumar2006}. Moreover, due to the limited
information content present in CW data, resolution is
entirely dependent on a tissue's optical properties and geometry; it
cannot otherwise be optimized \cite{Ntziachristos2005}.
These caveats can be overcome by employing frequency-domain or time- domain technologies.\\

\noindent {\bf Frequency-domain measurements}

Frequency-domain (FD) methods measure modulation amplitude attenuation and phase shift due to light propagation\cite{Nissila2002_RSI}. These FD measurements allow one to separate the contributions of absorption and scattering \cite{Intes2005a, Unlu2006} and image both fluorochromes concentration and lifetime bio-distributions \cite{Muraca2002_COCB, Godavarty2005_MP,
Bevilacqua2000_AO}. However, these advantages can be only achieved by multiple frequency systems \cite{Intes2005a, Bevilacqua2000_AO, Burcin2006}, where employing dense source-detectors arrays working at multiple frequencies in the range of a few MHz to 1GHz is extremely challenging \cite{Chance1998}.\\

\noindent {\bf Time-domain measurements}

Time-domain (TD) methods illuminate tissue with ultrafast photon pulses and resolve the arrival of the photons as a function of time at different locations around the tissue boundary. Photons collected shortly after the pulse experience significantly fewer scattering events and probe less volumes, while the late photons sample much greater volumes, which provide extra information useful for tomography. The temporal distribution of photons detected after propagation is known as the temporal point spread function (TPSF).
TD technology offers a high information content that is unmatched by CW and FD methods. The CW and FD data can be easily converted by integration and Fourier transform of the TPSF, where the latter has been applied to FMT due to the simplicity of the relationship between the fluorescence lifetime and the measured phase in FD data~\cite{Zhang2008, Soloviev2009}. Moreover, there are additional data generated from the TD measurements, such as the first second and third order moments that can provide enhanced separation between absorption and scattering \cite{Schweiger1999_PMB} or enhanced quantitative fluorescence tomography \cite{Lam2005}, which have no simple equivalent in the frequency domain.

%time gates..
The direct measurements of the time-dependent intensity using the gating technique is also of significant value. It is well established that resolution of the optical reconstructions can be improved by using early gates for both DOT \cite{Feng1993, Wu1997} and FMT~\cite{Niedre2008, Niedre2010} applications. Moreover, the time gated data is critical for multiplexing of fluorophores since it directly reflects the lifetimes on the measurements. Recently, studies have investigated the potential of performing FMT directly using time gates (TG)~\cite{Kumar2006}. In the case of lifetime multiplexed studies, unmixing algorithms applied to the decaying portion of the TPSFs have been employed to recover the fractional contribution of the different fluorescence components and perform FMT based on unmixed time-independent signals \cite{Kumar2008}. However, such methodology suffers from the same limitation as the continuous wave (CW) technique, i.e., limited resolution and requires unmixing algorithm where performances are not robust with low photon counts~\cite{Kollner1992}. Moreover, it is also limited to lifetime over 500ps \cite{Kumar2008}. Conversely, the use of discrete time gates spanning the full TPSF as the dataset for FMT should alleviate these drawbacks.

\subsection{\label{}Ill-posed inverse problem}
The ill-posed nature of the inverse problem leads to poor spatial resolution and reconstructions that are highly sensitive to noise. This ill-posedness can be alleviated by \emph{a priori} knowledge of spectral and structural properties of the imaged sample~\cite{Burcin2006}. Many strategies capitalizing on anatomical priors from a high resolution imaging modality such as MRI or X-ray have been devised for optical tomography, to improve the recovery of optical contrast~\cite{Azar2008}. Among these strategies, methods combining both anatomical and functional priors to guide the optical reconstruction offer superior quantitative performances for functional imaging applications~\cite{Intes2004}.

% MC
% parallel advancement
%Only recently the Monte Carlo method was accelerated thanks to the parallel computing advances. GPU and cluster use...
To summarize, an inverse algorithm based on a model that can accurately represent the light propagation in preclinical scenarios is in demand, with the potential of being incorporated with time-resolved data and various types of source-detector geometries. The Monte Carlo method, as the most accurate model in the time domain and in a wide optical property range, with outstanding flexibility for arbitrary geometries, can serve as an excellent candidate for this purpose. However, the computational burden of using the Monte Carlo method is further increased due to the photons required to achieve reliable statistics per time window. Until recently, time-resolved Monte Carlo approaches have focused on spectroscopic applications but not tomographic ones. Inversion strategies based on Monte Carlo models' CW components have been reported~\cite{Kumar2008a,%asymptotic cw
 Zhang2009}, and no formulation for direct time-resolved reconstruction based on Monte Carlo models has been yet reported.

\section{Thesis overview}
The principle objective of this thesis is to develop a highly accurate, efficient and flexible reconstruction approach for whole-body small animal optical tomography imaging. This approach will be based on the Monte Carlo method for accurately modeling light propagation in complex small animal models, thus achieving quantitative accuracy in reconstruction. It will be implemented in the time domain to benefit from the comprehensive information that the time-resolved data provide (simultaneously improving resolution and quantification in tomographic studies). Various excitation-detection strategies will be employed within this algorithm, utilizing the flexibility of the Monte Carlo approach. Furthermore, parallel computing methods will be incorporated to increase computational efficiency.

The integration of these developments will allow accurate and efficient reconstruction for quantitative 3D \emph{in vivo} imaging of functional and molecular information using temporal data and various source-detector settings, and fill an important unmet need in preclinical studies.

The remainder of this thesis will follow the outline below. In Chapter~\ref{chp:pmc}, the feasibility of performing time-resolved DOT based on the perturbation Monte Carlo method is investigated, for retrieving intrinsic optical property maps and functional parameter distribution with geometry guided source-detector settings. The validity of the model is experimentally established using an optical tomography system employing patterned illumination. Chapter \ref{chp:fluo} details the derivation of a computationally efficient model allowing for calculation of the weight functions for both absorption perturbation and fluorophore distribution while simulating only the propagation of excitation photons. Moreover, the gate information content associated with fluorescence lifetime is assessed to simultaneously reconstructing fluorophores with lifetime contrast. In Chapter~\ref{chp:comp}, a comparison of the perturbation Monte Carlo method to other Monte Carlo based reconstruction algorithms is presented in terms of time efficiency for time-gated tomographic applications, with an analysis of different parameters affecting the computational time. Lastly in Chapter \ref{chp:mesh}, a mesh-based Monte Carlo method for widefield illumination allowing for accurate boundary representation and fast calculation in tomographic applications is described. 